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Shape Sensitivity Analysis of Flutter Response of a Laminated Wing 

by 

Frederick D'Oench Bergen, Jr. 

Rakesh K. Kapania 
Aerospace and Ocean Engineering 
(ABSTRACT) 


A method is presented for calculating the shape sensitivity of a wing aeroelastic response 
with respect to changes in geometric shape. Yates' modified strip method is used in 
conjunction with Giles equivalent plate analysis to predict the flutter speed, frequency, 
and reduced frequency of the wing. Three methods are used to calculate the sensitivity 
of the eigenvalue. The first method is purely a finite difference calculation of the 
eigenvalue derivative directly from the solution of the flutter problem corresponding to 
the two different values of the shape parameters. The second method uses an analytic 
expression for the eigenvalue sensitivities of a general complex matrix, where the deriv- 
atives of the aerodynamic, mass, and stiffness matrices are computed using a finit® dif- 
ference approximation. The third method also uses an analytic expression for the 
eigenvalue sensitivities but the aerodynamic matrix is computed analytically. All three 
methods are found to be in good agreement with each other. The sensitivities of the 
eigenvalues were used to predict flutter speed, frequency , and reduced frequency. These 
approximations were found to be in good agreement with those obtained using a com- 
plete reanalysis. However, it is recommended that higher order terms be used in the 
calculations in order to assure greater accuracy. 
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1.0 Introduction. 


/./ [fUrodii^tion and Bac/iground 


Flutter, an aeroelastic instability, is a self-sustaining oscillation that involves the coupl- 
ing of inertial, elastic, and aerodynamic forces. Tail flutter was the earliest documented 
case of dynamic aeroelastic instability. It was discovered on the twin-engined Handley 
Page 0/400 Bomber at the beginning of World War I that the fuselage and tail had two 
principal modes of vibration. In the first mode, the elevators oscillated about their hinges 
ISO degrees out of phase, because the elevators were not attached by a torque tube. In 
the second made, the fuselage oscillated in torsion. The coupling between these two 
modes resulted in a self-sustained oscillation [lj. 

Modem aircrafts are subjected to numerous types of flutter phenomena. Classical flutter 
is associated with potential flow and involves coupling of two or more degrees of free- 
dom. Nonclassical flutter is concerned with oscillations due to boundary iaycr efTects, 
such as separcted flow and reattaching flows. Nonclassical flutter is more difficult to 
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analyze from a theoretical standpoint. Only classical flutter will be analyzed in this 
presentation. 

Flutter analysis capablities have been available for well over four decades. Loring [2] 
developed a general approach to the flutter problem in 1941. Yates [3) developed a 
modified strip analysis to analyze flutter characteristics for finite-span swept and un- 
swept wings at subsonic and supersonic speeds in 1958. This method is still used today 
to calculate the lift and moment forces. For example, Landsberger and Dugundji [4] used 
these expressions, with a modification for camber effects given by Spielberg [5], to study 
the flutter and divergence of a composite plate. The present day computers allow still 
more complex aerodynamic programs to be developed and used. To name a few, the 
unsteady vortex-lattice method developed by Strganac and Mook [6], the double lattice 
method developed by Waldman (7], and the MCAERO code developed by Hawk and 
Bristow [81. 


The dynamics of any physical system is important from a designer's view point. The 
onset of flutter in a modem aircraft will involve large oscillatory distortions of the 
str uctural components. Therefore, it would be advantageous to the designer to have a 
tool which can be used to predict the changes in flutter with the changes in basic shape 
parameters. 

Sensitivity analysis was first recognized as a useful tool for assessing the effects of 
changing parameters in mathmatical models of control systems. The gradient based 
mathematical programming method used in optimal control and structural optimization 
furthered the development of sensitivity derivatives, because sensitivity derivatives are 
used in search directions to find optimum solution [9,10,1 1,12]. Inefticicnt optimization 
programs lead to further interest in efficient calculation of sensitivity derivatives. As a 
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result, sensitivity analysis has become a versat.le design tool, rather than just an instru- 
ment of optimization programs [13j. 

Shape sensitivity analysis of any physical system under aeroelastic loads can be impor- 
tant from different points of view: 1) to understand and predict the system's response 
and 2) to optimize the response of the system for a set of physical constraints. In order 
to predict or optimize the response the sensitivity derivatives must be calculated. This 
can be done by a finite difference calculation or analytical calculation of the sensitivity 
derivatives. Analytical sensitivity analysis has found increased interest in engineering 
design. The analytical derivatives eliminate the uncertainty in the choice of step size, 

which if too large can lead to truncation errors and if too small can lead to ill- 
conditioning. 

Adelman and Haftka [13] have shown that structural sensitivity analysis has been avail- 
able for over two decades. Structural sensitivity analysis has been sufficient in the past 
because sizing variables such as plate thickness and cross-sectional areas effect the mass 
and stiffness properties of the airframe but, not its basic geometry. Therefore, aero- 
dynamic sensitivity analysis capability has been limited in development until recently. 

For example, Rudisill and Bhatia [14] developed expressions for the analytical deriva- 
tives of the eigenvalues, reduced frequency and flutter speed with respect to structural 
parameters for use in minimizing the total mass. However, this method is limited because 

the structural parameters are sizing vanables such as cross-sectional areas, plate thick- 
ness and diameters of spars. 

Pedersen and Seyranian [15], examined the change in flutter load as a function of change 
in stiffness, mass, boundary conditions or load distribution. They showed how sensitivity 
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analysis can be performed without any new eigenvalue analysis. The solution to the main 
and an adjoint problem provide all the necessary information for evaluating sensitivities. 
Their paper mainly focuses on column and beam critical load distributions. 

Hawk and Bristow [8] developed aerodynamic sensitivity analysis capabilities in suberit- 
ical compressible flow. They first analyzed a baseline configuration, and then calculated 
a matrix containing the partial derivatives of the potential at each control point with 
respect to each known geometric parameter by applying a first-order expansion to the 
baseline configuration. The matrix of partial derivatives is used in each iteration cycle 
to analyze the perturbed geometry. However, this analysis only handles chordwise per- 
turbation distributions, such as changes in camber, thickness and twist. A new approach, 
which is still under development, has been proposed by Yates [16] that considers general 
geometric variations, including planform, for subsonic, sonic and supersonic unsteady, 
nonplanar lifting-surface theory. 


Barthelemy and Bergen [17] explored the analytical shape sensitivty derivatives of a wing 
static aeroelastic response. They found that the second derivatives of the wing's 
aeroelastic characteristics, such as section lift, angle of attack , rolling moment, induced 
drag, and divergence dynamic pressure, for subsonic subcritical flow, with respect to 
geometric parameters are small enough to be well approximated by sensitivity based 
linear approximations. These approximations are valid within a range that are useful to 
the designer in the initial design phase. 

This work is an extension of previous work preformed by Barthelemy and Bergen [17], 
and details the theoretical and computational derivation of a method to obtain the sen- 
sitivity of a wing flutter response to changes in its geometry. Specifically, the objective 
is to determine the derivatives of flutter speed and frequency with respect to wing area, 
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aspect ratio, taper ratio, and sweep angle. Chapter 2 gives a brief description of the 
structural formulation whtch was onginally formulated at the NASA I.angley Reasearch 
Center by dies ( 18 |. The program is based upon a Raleigh- R,.z formulation in wh.ch 


the displacement functions are made 
formulation is presented in Chapter 3. 


up of polynomial expressions. The aerodynamic 
The expression for lift and moment are derived 


from potential now theory and have been modified to account for finite span. In 
Chapter 4, the formulation is validated using examples found other works. Once there 
is sufficient confrdence in the flutter speed prediction capabli.ies, is shown that the 
sensitivity analysis will predict the flutter speed for changes in the structural parameters. 
Three different approaches are used to obtain the sensittvity of the flutter speed, fre- 
quency, and reduced frequency to various shape parameters. These include: 


(i) a purely numerical approach using the finite difference method; 
(it) a semt-analytical method that uses an analytic exf.ession given 
by Murthy and Haftka (19|, for calculating the sensitivity of the 
eigenvalues of a generalized eigenvalue problem and a 
finite difference approximation of the derivatives of the 
aerodynamic, mass and stiiTness matrices with respect to the 
geometric parameters; and 

(iii) an anal >' tic a PProach that uses the analytic expression for 
calculating the sensitivity of the eigenvalues and 
an analytically derived expression for the sensitivity of the 
aerodynamic matrix with respect to geometric parameters. 

Chapter 6 gives conclusions and suggestions for future work. 
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2,0 Structural Model 


2,7 Introduction 


This chapter presents the structural formulation which was originally programmed by 
Giles (18]. The program is based upon a Ritz solution technique using the the energy 
functionals for a laminated plate which includes the bending and stretching of the refer- 
ence surface. This program is capable of analyzing unsymmetric wing box sections aris- 
ing from airfoil camber, laminate sequences, or different thicknesses in upper and lower 
cover skins. Thermal loadings can be represented as temperature distnbutions over the 
planform of the cover skins. 

The aerodynamic formulation restricts the chordwise length to remain straight during 
oscillations, therefore only high aspect ratio wings will be analyzed in this presentation. 
Only bending and torsional deformations are considered in this analysis procedure. In 
the theory of plates, the Kirchhoff assumption is made that lines normal to the reference 
plane remain straight and normal after deformations. 
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The planform geometry of the wing is represented by a trapezoidal segment. To repre- 
sent cranked wing boxes, multiple trapezoidal segments can be defined. Each segment 
has a separate local coordinate system. The local coordinates are nondimensional such 
that <* refers to a fraction of the chord and ij a fraction of the span as shown in Figure 
1. 

Figure 2 illustrates a possible geometry of the wing box and coordinate system. The 
mid-camber surface is measured from a reference plane. The distance, Z c , is represented 
as a polynomial in the global coordinates x andy. 

Z c (x y) = z 00 4- z 10 jt 4- z 20 x 2 + z 0l y 4- ... 4- z mrr x m y n [2.1.1] 

The depth of the wing box, H(x,y) is measured from the midchord surface and the 
thickness t(x,y) of each laminate is also defined by polynomial expressions x and y. The 
expressions for the depth and thickness are: 

H(xy) = H qq + H l0 x + H 2 Q x 2 + H 0l y + ... + H mr x m y n [2.1.2] 

t(xy) = r 00 4- / 10 Jf + t 2Q x 2 4- roiy + ... 4- t mn x m y n [2.1.3] 

In the present formulation, the depth of the wing box and thickness of the skins is as- 
sumed to be constant throughout. The coefficients of equations (2.1.2) and (2.1.3) are: 


Hqq =» .100 meter * 

r 0 o = .002 

meters 

Hmn ” 0 0 fmn-0.0 

fit ® 1 5 

n * 1 


The wing box used for this presentation is shown in Figure 3. The box section chordline 
is straight, therefore Z e ( x,y) is zero for this presentation. 
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2.2 Displacement Function 


The Rayleigh- Ritz formulation assumes a deflection shape for the wing structure. This 
deflection shape is a linear combination of n assummed displacement functions. The as- 
sumed displacement functions are specified as products of polynomial x-direction and 
y-direction global coordinates. The deflection equation can be written as: 




n=0 m=s2 


^max 


\m 


[ 2 . 2 . 1 ] 


where N is restricted to five terms and M is restricted to six terms in order to prevent 
numerical difficulties in the manipulation of the matrices. W{x ,y) is the transverse de- 
flection of the reference line , and C„, the set of unknown coefficients. 

The deflection equation can also be written as: 

np 

W ( x >y)~ X y «' (jf C « [2-2.2] 

<=l 

where y, ( x, y) is the nondimensional displacement function and np is N times M. All 
the assumed displacement functions satisfy the geometric boundary conditions for a 
cantilever plate. The constants x^ and^ are defined in magnitude as the root chord 
and span respectively. The displacement functions are nondimensional quantities in or- 
der to prevent numerical difficulties in manipulation of the matrices. 
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2.3 Energy Expressions 


The strain energy U for a symmetric anisotropic laminate is given as: 

+ v„ + d 22 iv 2 „ 

/t 

+ 4 Dlt"' J ,»' Ja + 4D a ir j , f iy„ + 4D u tvi, )JA [2,3.0 

Where a comma denotes panial dt.Teremiation and D, t are the bending stifTncss terms. 
The kinetic energy T for the plate is: 


T 



mW 2 dA 


[2-3.2] 


where m is the distributed mass per unit area IV is the transverse velocity at or andy, and 
A is the area of the plate. 


The Ritz soiuuon procedure is used to determine the numerical values of the set of un- 
known coefficients, C,'s. To use the Ritz solution procedure, and to integrate these 
equations the global x and y were transformed into a local system. This was achieved 
by the following transformation. 


x - e + a’ ‘ + (f- e)rj + (c - 


[2.3.3] 


y =* g‘ + Brj 


[2.3.4] 
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The coefficients a', B, c, e,f, g' are dimensions of the local coordinate system as shown 
in Figure 1. 

The Ritz solution procedure is based upon the second variational priniciple. The second 
variational priniciple states that the total energy of the structure can be written in terms 
of any kinematically admissible function. If the total energy, E{y), in terms of the 
kinematically admissible function y, is stationary at y with E(y) - 0, then the func- 
tion represents the real modes of the structure. If the function represents the real modes 
of the structure then it can be used as a valid approximation to the exact solution. The 
total energy of the structure can be written as: 

E “ u ~ T [2.3.5] 

The extremum principle states that the energy is stationary with respect to Cj i.e. 

dc" = ® P r °duces a system of n simultaneous equations. In matrix form then the 
equations are expressed as: 

[[AT) - o) 2 [A/]]{C) - {0} [2.3.6] 

where [KJ and [M! are the stifTness and mass matrices resulting from the partial differ- 
entiation of the potential and kinetic energies, recpectively, with respect to the general- 
ized coordinates. 
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Acyod\jtanii£ Forces 


^ chapter presents the formulation of the aerodynatmc coemcien. ms<nx. An 

“ C0 ~- 2 -“* «* was calculate the aero- 

Onamrc coeE, cents. Thts method was fust developed by Bam.by, et. al.. |20| and was 

modified by rates to tndude the effects of fnnte span. Lift and moment forces are 

' mcd a.ong the midchord and acting upon sections perpendicular to th.s midchord line 
(called the reference line hereafter). 


The now field is represented by a uniform stream (non-cmculatory component, super- 
Posed by the dtsturbance-ve.ocy d.stnbutton (circulatory component, whtch wr„ mode, 

,h ' ° fthe P0!U,0n “ d 0fth ' -** " - "-deled such that the ennditton 

Of tangential How at the mng surface is met. 


The position of each point of the 


segment taken perpendicular to the reference line is: 
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Z~h + x'd 


C3 .1.1] 


Where h is the deflection of the midchord, positive up, 9 is the rotation of a section 
measured perpendicular to the midchord positive leading edge down, and x' is the per- 
pendicular distance from the reference line, measured positive back. 

The lift and moment forces are written in terms of a circulatory and non- circulatory 
components. As taken from Yates' modified strip analysis (3), the lift force (positive up) 
per unit length of the wing, is given as: 

L = - npb 2 [h + V n 9 + V n o tan A - ba(9 + V n r tan A)]„ c 

-{Ci aM pV n bC{k n )Q d ] c [3.1.2] 

where subscript nc indicates the noncirculator)' component and c indicates the 
circulatory component. The air density is p, b is the halfchord, V n is the air speed meas- 
ured perpendicular to the midchord. C,^ H is the section lift coefTiuent which is constant 
along the span, A is the midchord sweep angle, and a is the off-set of the reference line 
from the midchord line, assumed to be zero in this presentation. The local bending slope 
of the reference axis is represented by a and the local rate of change of twist is repres- 
ented by r. C(k„) is Thcoderson's circulation function. The reduced frequency, k m is a 
measure of the amount of circulation in the flow field and is used in the circulation 
function and is given as (10): 


0 165 0.335 


0.0455 . 0.3 . 

~~~k — ' 1 “T"' 


[3.1.3] 
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Similarly, the expression for the moment M, per unit length about the midchord (posi- 
tive leading edge down) is : 

A/ = -npb*[{ -5- + a 2 )(0 + V n x tan A) + npb 2 V n (h + V n o tan A) + npb } a{fi + V n a tan A) 

O 

4- np y 2 b 2 {d - abx tan A)]„ c 

[ -2*i> v„b \ y - (« - i„)C(*) -%-]&], D-i-fl 

A A. 71 

where a cn is the distance to the aerodynamic center from the reference line, and Q d is the 
downwash distribution defined by: 

C, 

Qj = h + V n d + V n a tan A 4- b{ ~~ + a cn - a)(d -r V n x tan A) [3.1.5] 

The airfoil geometry with typical dimensions is shown in Figure 4. 
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The lift and moment forces can. be reduced to the following expressions assuming the 
wing is undergoing infinitesimal harmonic oscillations about its steady-state position 
[ 20 ]. 



L~ — npb*o> 2 {B cft h + B cd 9) 

[3.1.5] 


M = - npb*(o 2 {B ah h + B aQ d) 

[3.1.6] 

The coefficients B cH B c9 , B A „ and B^ are defined as: 



B ch = \ A ch~ A ch tan A 

[3.1.7] 


S c6 = A c* + Y bA ci tan A 

[3.1.8] 


B ah = \ A ah~ -JJ A ah tan A 

[3.1.9] 


B a0 = A aa + ^ A ar tan A 

[3.1.10] 

and; 

4 t , , 

^"' ,+ 

[3.1.11] 


. , . 2<-a*JC ( „ „ _ 2 <WC,„ 

A c* a + k + •. . (1 2a) + 

*« *;2* 

[3.1.12] 


, /a . 2C(A rt )C /a ^ / 

^ CT *„ + k n 2n ( a « d) 

[3.1.13] 
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np 

°( x •?)- ,y) cos A sin A 

<=l 

+ Yi'iyi* <y){ cos 2 A - sin 2 A) - y w (x ,y) cos A sin A) [3.1.19] 

Here the superscripts - represent the fact that x, y are not any arbitrary values of x and 
y , but x , y arc the coordinates of the reference line in x and^y. Futhermore: 

x =y tan A [3.1.20] 

Therefore h y 0, t, and c are functions of y only, where y( = — — — ) is the distance of 

cos A 

the point x f y from the origin along the reference line. For simplicity, A, d , t, and a can 
be written as: 


np 

K y) - X c < Y ( y> 


1=1 


np 

l=l 


np 


- £c, r Hy) 


1=1 


np 

tOT-Vc, 7,,-OT 


1=1 


[3.1.21] 


[3.1.22] 


[3.1.23] 


[3.1.24] 
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The lift and moment forces can then be written in terms of the displacement functions 
y,(y) and the unknown coefficients, C r Based upon equations 3.1.5, 3.1.6, and 
3.1.17-3.1.24, the lift and moment can be written as: 

np np 

L = — npb U) VmCv) On ~ ^ Vm, y^) On) ^ ch 

m~ 1 n m = 1 


np np 

+ A cx ^ y m% j(y) C m +bA c < ! tan A ^ y m§ jjtIk) [3. 1.25] 

m=l 1 


np np 

M = — Trpfr { ~ V m (> , Qn ~~ ^an ^ jr(V ) 

m=l n «=1 

np np 

+ ^«X^^ Cm+WflT tanA j]y w ,53?(nC , „ [3.1.26] 

m=l "»=1 

Having obtained the lift and moment forces, in terms of the displacement function and 
unknown coefficients, attention must be focused on the principle of virtual work since 
these forces are non-conservative. At this point forward, y(y) will be represented as y for 
ease in derivations. 
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3.2 Virtual Work 


The lift and moment forces are non-ccnscrvative forces, therefore the principle of virtual 


work must be employed. The definition of virtual work gives: 


1 Wnc = ( ‘tehdy - J C ; 


[3.2.1] 


7=1 


where <5/i, and -50 arc the virtual displacements and Q j arc the generalized forces. The 
virtual displacements are defined as: 



[3.2.2] 


np 


se-Yy^Cj 


7=1 


[3.2.3] 


The virtual work can be expressed in terms of generalized forces and displacements as: 


np 


[3.2.4] 


This can be written as: 
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• r~ 


np 


np 


SIP 


nc 


Vj)C t - y- tan A(7 ( - y)C 

;=i u * , * 

i = l 


+ (?f, x V,) C,rl ca + b tan A (y u - y y )C,v} CT dy 5 C 




- J up* [ j- (v, jj. ;) C, - f- tan Aft. ? y, r ) CJ, 


1=1 


+ yi. I y>. x + * tan A(y (> xJT Vj, x)^i A an ] 4F«>C. 


The generalized forces are defined as: 


np 


v c, 


1=1 


The aerodynamic matrix then becomes: 


tan Aft.,- y ; )M rt 

+ y«. f yj A c* + * tan Afy (> ^ 7 y M CT } ^ 



^P^ 4 {[ J (y, Yj, j) - y- tan A(y, <jr y,- r )K„ 
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+ Vi, x Vj, x'^aa. ^ xy Vj, 


[3.2.7] 





TOT- 


w 



1 

' 


The integral equations for the aerodynamic matrix are obtained by substituting the ex- 
pressions for the aerodynamic coefficients, and the displacement functions y 's into the 
expression for the virtual work of the non-conservative forces. The integration is per- 
formed numerically by a 15 point gaussian quadrature numerical integration scheme [23]. 
The roots of the Legendre Polynomial and the weight factors used in the numerical in- 
tegration are given in Table 1. 
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Tabic 1: Abscissas and Weights in Gaussian Quadrature, n* 15 


Abscissas 


0.000000000000000 
+/-0.201 194093997435 
+ /-0. 394 1 5 1 347077563 
+/-0.570972 172608539 
+ /-0.7244 1773 1360 170 
+/-0.848206583410427 
+/-0.937273392400706 
+ /-0.9879925 18020485 


Weights 

0.202578241925561 

0.198431485327111 

0.186161000115562 

0.166269205816994 

0.139570677926154 

0.107159220467172 

0.070366047488108 

0.030753241996117 




3.3 Flutter Analysis 


The V-g method [24] was used in computing the flutter speeds. This method introduces 
a structural damping coefficient g into the equations of motion. Neutral stability (flutter) 
is attained for a given velocity, when the damping of the structure goes to zero. As- 
suming harmonic motion the equations of motion become: 

[CW + ig) ~ « 2 [iV]]{C} = c o\a] {C} [3.3.1] 

In the absence of non-acrodynamic forces, the resulting generalized eigenvalue problem 
can be written as: 


CCW -[*]]{$}«{<>} 


[3.3.2] 


B is a generalized complex matrix, X is the eigenvalue and {C} is the eigenvector. The 
eigenvalue /I is defined as: 


1 + ig 


The generalized complex bar matrix bar [B] is defined as: 

[*] - [*r‘[.Vf + A] 


[3.3.3] 


[3.3.4] 


The flutter speed perpendicular to the midchord is defined as: 


i . (ijb . 

V n-~T cos A 


[3.3.5] 
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Where b is defined as the halfchord perpendicular to the rrudchord reference line. The 
flutter frequency bar is w, and k„ is the reduced frequency. 
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iQ Validation of Analysis Program 


4. / Introduction 


This chapter gives comparisons between results found in literature and the present code. 
The validation process is broken down into four sections. To validate the stiffness ma- 
trix, the static deflections are checked. The mass matrix is verified by comparing the 
natural frequencies of vibration with isotropic as well as composite materials. Once suf- 
ficient verification of the structural model is complete, the static aerodynamic loads are 
checked for divergence of swept and unswept wings. The dynamic aerodynamic matrix 
is verified by comparing the flutter frequencies and speeds with results found by Castel 
and Kapania (21) who developed a simple element for the aeroelastic analysis of lami- 
nated wings. Their formulation allows for unsymmetric laminations, arbitrary geometry 
including chord and thickness taper, and multiple sweep angles. 
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The comparisons made in this presentation are for a wing consisting of top and bottom 
flat laminated skins rigidly connected as shown in Figure 3. For all isotropic compar- 
isons, the material properties of the rectangular box beam were taken to be Aluminium: 


£ n = £ 22 , 6.89 48 x-lO 10 ^- v,2 » 0.30 

rn 

The material properties of the material used in the laminated wing are: 

£ n -6.9* 10 10 -^-, ^22 * 5.0 x 10 10 v 12 = 0.2S 

m m 

Gu-UxlO 10 -^ p mal = 2.71 * lu 3 — y 
m m 

It should be noted that the material properties were arbitrarily chosen. 


4.2 Static Deflections 

The static deflections are compared with Beer and Johnston [22], and Castel and 
Kapania [21]. For tip loading of a uniform cantilever beam, Beer and Johnston (22] give 
the tip deflection to be: 


Z- 


uL 

2 EI 


(positive up) 


[4.2.1] 


Comparisons with the results obtained using the EPFAC (Equivalent Plate Rutter 
Analysis Code) are given in Table 2, for different aspect ratios, AR. 
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Next, the tip defection under the same loading condition was tested against results of 
Castel and Kapania [22] for the wing in Figure 3. In this formulation, Castel and 
Kapania used four finite elements to model the structure. Each beam element has 
twenty-four degrees of freedom, twelve at each end. Results for zero sweep and variation 
in fiber angle are given in Figure 5. Excellent agreement between the finite element code 
and EPFAC exsists for the entire range of values of the fiber angle. 

Figure 6 gives a clear comparison between the two different models. EPFAC is termed 
a skewed wing while the model of reference 21 is a rotated wing. The differences in de- 
flections obtained using the different modeling of the wing structure can be accounted 
for as follows. The rotated wing's boundary conditions are such that the analysis pro- 
cedure assumes an eflective root and tip perpendicular to the reference line, regardless 
of sweep angle. As a result, the static deflections arc independent of the sweep. The 
skewed wing assumes the root and tip parallel to the flow. As the wing is swept forward 
or backward, the torsional and bending stiffnesses become coupled. Static deflection 
comparisons for different sweep angles of the composite material are given in Figure 7. 


4.3 Fi ‘ee Vibration Analysis 


For the wing under consideration, the natural frequencies obtained from the present 
code EPFAC was compared with those obtained using the expressions for a uniform 
Bernoulli-Euler Cantilever Beam [25]. The comparisons of these values are shown in 
Table 3. The results show excellent agreement for different aspect ratios. 
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EPFVvC was compared with Ref. [21J, for an unswept wing having different ply angles. 
The results are plotted in Figure 8. The natural frequencies of forward and aft swept 
wings were also tested. The results are shown in Figure 9. The differences between the 
two sets of results at larger sweep angles are due to the differences in modeling of the 
boundary conditions, in the two studies as discussed previously in section 4.2. 
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Table 2 : Comparisons of Tip Deflections (meters) for Different Aspect Ratios 


AR 

EPFAC 

Eqn. 4.2.1 

%diff 

5 

- 0.04748 

-0.04646 

2.14 

10 

-0.18995 

-0.18587 

2.15 

20 

-0.75968 

-0.74342 

2.14 

30 

-1.71032 

-1.67408 

2.12 
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SWEEP ANGLE (DEG) 


Fl|urt 7 Tip Deflection Comparison for Different Sweep Angles 
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FIBER ANGLE , DEG 


Figure 8 Variation of the Natural Frequencies with Fiber Angle 
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Figure* 


Variation of the Natural Frequencies with Sweep Angle 








4.4 Divergence and Flutter 


Having gained confidence in the structural model, the results for aeroelastic response of 
swept and unswept wings were compared with the results from two different codes: ( 1 ) a 
code written by Barthelemy and Bergen [17), and (ii) a code written by Castel and 
Kapania [21]. Barthelemy and Bergen used Wcissinger's L-Method to obtain the static 
aerodynamic loading matrix. 

In Wcissinger's L-mct od, the flow around the wing is modeled by a lifting line of 
vortices bound at the wing quarter-chord line. A no-pentration boundary condition is 
specified at na control stations along the win*’ The placement of the control points de- 
termines the spanwise distribution of ' . modification to the position of the down- 

wash distribution was employed to account for the slope of the C, vs a curve less than 
2rr. 

The wing dimensions used for comparison with the present method are: 

S- 20.0 m 2 /f/?=10.0 tp- 1.0 

where tp is the taper ratio of the half wing. The results are shown in Figure 10. An 
excellent agreement is observed between the two sets of results. 

Castel and Kapania [21] used Yates' modified strip method [3] to obtain both the static 
and dynamic aerodynamic loadings for swept and unswept composite wings. For dif- 
ferent values of sweep the fiber orientation was varied and the results for the divergence 
and flutter speeds were obtained with EPFAC. Figures 11-15, show the comparison 
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between the results obtained using the EPFAC code and that of Castel and Kapania for 
divergence and flutter speeds for different va’.jes of sweep and fiber orientation. The 
differences at large sweep angles arc due to the different models used in the two studies. 
As discussed previously in the section on deflection in this chapter, the model used in 
Reference (21] is a rotated model, compared to a skewed model used in EPFAC 
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SWEEP ANGLE (DEG) 


Figure 10 Variation of Divergence Speed with Sweep Angle 
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Figurt 12 Divergence and Flutter Speeds of a 30 * Sweptforward Box Beam 
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5.0 Sensitivity Analysis 


5 . 1 Introduction 


This chapter presents the calculation for the sensitivity of the flutter speed, flutter fre- 
quency, and reduced frequency to geometric shape parameters namely: (i) aspect ratio, 
(ii) surface area, (iii) taper ratio, and (iv) sweep angle. The sensitivity calculations require 
the sensitivity of the aerodynamic, mass and stiffness matrices with respect to various 
shape parameters. A key objective of this study is to check the viability of calculating 
the desired derivatives using an analytical approach. It was decided to analytically ob- 
tain the sensitivities of the aerodynamic matrix. The analytical derivatives eliminate the 
uncertainty in the choice of step size which if too large can can lead to truncation errors 
and if too small can lead to ill-conditioning. 

The derivation of the analytic derivative of the aerodynamic matrix, is given in this 
chapter. To validate the expressions for the eigenvalue derivatives, these derivatives are 
calculated using three different methods. The first method is purely a numerical ap- 
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proach that uses a finite difference approximation to find the eigenvalue derivatives. The 
second method is a semi-analytic approach, because the derivatives of the aerodynamic 
matrix are found using finite difference approximations, and then using the expression 
for the derivative of the eigenvalue as given by Murthy and Haftka [19]. The third "an- 
alytic" method uses an analytically derived derivative of the aerodynamic matrix, along 
with the eigenvalue derivative expression given by Murthy and Haftka [ ' 9]. 


5.2 Eigenvalue Derivatives and Solution Procedures 


In the first method of calculating the derivatives, the flutter problem is solved twice and 
the derivatives of the eigenvalues are approximated using a forward finite difference 
scheme. 

The second and third methods use the expression for the eigenvalue derivative as given 
by Murthy and Haftka [19]. The expression is derived using the main and the adjoint 
problem. The main eigenvalue problem is: 

C[^-[5]]{e r } = {0} [5.2.1] 

where e, is the right hand eigenvector, and (B) is defined in equation 3.3.4. Similiarly the 
adjoint problem is: 




[5.2.2] 
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where e t is the left hand eigenvector. The eigenvalues in both cases are the same. Taking 
the derivative of equation 5.2.1 with respect to a shape parameter p, yields: 


[[/] 



£[g] 

3p s 




de r 

5p s 


= 0 


[5.2.3] 


Multiplying the above equation by left eigenvector gives: 


M'Mf- 


3[£L 

dp s 


1M + {‘/J'CCW 


[*]] { ■— } = 0 

c Ps 


[5.2.4] 


Equation 5.2.2 reduces this to: 


{*/}'[[/] -If- 

dp s 


£[ 5 ] 

dp s 


]W=0 


[5.2.5] 


For the eigenvalue, the eigenvalue derivative is expressed as: 


8X l 

8 Ps 


, i x t dlB] , 
{ *' } (<? - } 


[5.2.6] 


Recalling equation 3.3.3, the eigenvalue derivative in terms of flutter frequency, damp- 


ing, and there derivatives is written as: 


dX 

d Ps 



3 

CJ 


+ 


r 2 - 8 co , 

[~ — CO -2(0- — g] 
cp s op t 


4 

0) 


[5.2.7] 
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To obtain the derivative of the generalized complex matrix [B], the derivatives of the 
aerodynamic matrix [A]; of the inverse of the stiffness matrix [AT} -1 , and of the mass 
matrix (M) are needed. 


£ 5 ] 

d Ps 


w 1 

dp s 


cc-w] + ixn + ixr'c 


cp s 


Mi 

dp 


[5.2.8] 


The derivatives of the mass matrix, and the inverse of the stiffness matrix are obtained 
using the forward finite difference method. A study was first conducted to obtain an 
appropriate step size for the finite difference calculations. The results of this study are 
shown in Table 4. These step sizes indicate that the calculated derivatives are stable as 
indicated in Table 4. 


The derivatives of the aerodynamic matrix [Aj with respect to a geometric shape pa- 
rameter are obtained using two different methods: (i) finite difference method; and (ii) 
analytic method. The calculation of the sensitivity of the aerodynamic matrix [A] is made 
difficult by the fact that this matrix depends upon shape parameters p, and also on the 
reduced frequency k n . The reduced frequency is not really an independent variable, as its 
value for a new value of p ^ = tf‘ d + Ap,) should be such that the imaginary part of the 
eigenvalue corresponding to the perturbed configuration should be zero. 

In the finite difference method, the sensitivity of the aerodynamic matrix [A] is obtained 
as follows: 


d[A] _ IA(P, + Aft. + &*„)] - [% k n y] 
d P, " *P S 
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To obtain the value of &k„ an iterative procedure, as described in Figure 16, is used. 
As a first step, Ak n is put equal to zero, and the sensitivity of the eigenvalue is obtained. 
Obviously, the imaginary part of the new eigenvalue thus obtained will not be zero. This 
fact is used to obtain the value of Ak n as explained in the following. 

The change in the damping coefficient g, a function of both the shape parameters and 
reduced frequency, can be written as: 

dg= if; dp s + jr n dk » f> 2 - 10 ] 


At flutter speed: 

g - 0 , dg = 0 


Therefore 


Jg_ 

dk n dp s 

[5 ' 21| i 

a*. 

C 

Note that — was already obtained during the calculations of the flutter speed. The 
values of reduced frequency are varied in the mtial problem to compute the value at the 
point the damping coefficient goes to zero. Therefore, ~~ is easily computed by a for- 
ward finite difference scheme. The value of-0- is obtained from the imaginary part of 
the sensitivity of the eigenvalue obtained in the first step. This is computed directly from 
the eigenvalue derivative. Recalling equation 5.2.7: 
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II 

5p s 


= imag ( 


d). 

dp s 


) w 2 + 



g 


cj 


where co is obtained from the orginal flutter problem and; 


doj 

5 Ps 


„ dk . 3 

real{ — — ) oj 

d Ps 

2 


[ 5 . 2 . 13 ] 


eg 

If ‘ s not w ithin a tolerance of 10' 4 , the aerodynamic matrix [A] is recalculated at a 

new value of the reduced frequency while also keeping the same perturbation in the 

dk n 

shape paramete-. Once -^j-is known, an approximation to the value of k n , correspond- 
ing to flutter speed for a new value of p t , is obtained using a Taylor's series expansion: 


b new _ i,old 

*n — K n 


4- 


dK 

d P: 


A Ps 


[5.2.14] 


This process is repeated until the tolerance is met. 


In the 'analytical' method, the sensitivity of the generalized comp. ex matrix [B] is ob- 
tained in a sirruliar fashion except the sensitivity of the aerodynamic matrix [A] is com- 
puted analytically. This is expressed as follows: 


4^] _ d[A] + d[Al dk^ 

d Ps " d Ps dk n d Ps 


[5.2.15] 


Both ^ - and — rr - arc derived in section 5.4. The value of — — - is computed from 
op, dp, dk v 

the eigenvalue derivative, as explained above. In the first iteration,—— is assumed to 

a P, 

be zero. The aerodynamic matrix derivative is computed and combined with the deriva- 
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tives of the mass [M], inverse stiffness matrices [K] -1 to form — The eigenvalue 


derivative is then computed. The value of — — is calculated as described above. Figure 

dp, 

17 shows a detailed flow chart of the analytic method. 


The derivative of the flutter speed is: 


evj_ 

dp s 




CD 


dk n 

] -cut-—-] 

cp s dp s 


[5.2.16] 


The flutter speed at the new configuration is then computed using a truncated taylor 
series. As before: 

vf-vf* [ 5 . 2 . 17 ] 

To calculate the analytic derivative of the aerodynamic matrix, the geometric parameters 
must first be defined. 

i 

I 

f 

5.3 Geometric Parameters j 

t 

The Aerodynamic matrix is composed of geometric parameters which are a function of 
surface area, aspect ratio, taper ratio and sweep angle. These parameters are the span, 
root chord, tip chord, halfchord, the global x and y coordinates, the displacement func- j 

i 

tions, and their derivatives with respect to global x and y coordinates. These are ex- f 

pressed as: 
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B = Jar s 


[5.3.1] 


C 


25 


(1 + tp)sjAR S 


C, 


2 Srp 


(I +tp)jARS 


crd(i)^C r + {C t -C r )^~l 


t= 1 15 


•»(') = -—+></) tan A /=1 15 


A‘) - 


Bri(i) + B 
2 


‘= 1 ..... 1 5 




crdii) 

2 


cos A 


1 15 


[5.3.2] 

[5.3.3] 

[5.3.4] 

[5.3.5] 

[5.3.6] 

[5.3.7] 


Where b is the halfchord measured perpendicular to the reference line , AR is the aspect 

ratio, 5 is the surface area. C, is the root chord, C, is the tip chord, crd(i) is the chord 
length at station i , and tp is the taper ratio. 


5.4 Analytic Derivatives 


For general purposes, p, will represent the shape parameters, where s - I is the surface 
area. , - 2 is the aspect ratio, s - 3 ,s the taper ratio, and s - a is the sweep angle. The 
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Aerodynamic matrix, given in equation 3.2.7, must be differentiated with respect to the 
shape parameters and the reduced frequency. In order to avoid the use of Lcbnitz's Rule, 
and facilitating numerical integration using Gauss integration, the intcrgrals must be 
changed over to local coordinates. Equation 3.2.7 is written as. 


A ji = - 1 npb^ai 1 [[ y (y, y y ) - -j- tan A (y,j Vj)]A ch 


B 


+ x Yj^cct + b xy 2 COS A ^ 


-r 

j -i 


npb (o [[ ^ (>'j '/j, x) k tan Mv,j y ;i i)1 A ah 


[5.4.1] 


B 


+ A. x Yj. x A aa ^ tan A (y (< xy y j^)A ax ] ^ drj 


g 

where dy is replaced by — — drj. 

2 cos A 


The derivative of the aerodynamic matrix with respect to a general shape parameter, 
D , B 

representing B = — — — . is: 


AA /< lr f 1 rnk c'b .. .2 ^1 - .2- . J 

— U2 b — h>/ + b —y,* b A c 


fib- db . 


ib' r Kf - 


ib 1 . l h 


r J i U {Ju - - 4 i U '% S — * , 1 Lf — J . A 

K CPs ' 1 '*'' A+ k n dp s /j tan A + k n y ‘-> dp, 


ib' 


d A 3 - - 


2k. 


- iU - - . _ _ , T -I , 

y '-pj co. ,2 y bpj tan A ] A ch 


k n cos\\ ' c Pt k n 
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+ (4b 


. /-m.2 db - - 3 **.* _ 3_ < 5 / 

+ (3 ‘ » —yj + >> a c . 


i db - 
$Ps 


4 *7 


Svi 


y,. ~i t tan A + b -f— y, tan A + tan A + 7, -y, -fi 


cp s 


rr, 2 \ Y ‘' XyYj 3 d ^ A °~ 

COS /V L Ps 


, n ,2 db __ 3 _ 1 ^.7 

+ (3J ^’W 7 .f+ > — /,.? + * >,-£-) 4 




GPj 


ah 


r 4ib^_db- _ . /ft 4 ^Vt.y ;a 4 cy - 

C *. ,anA+ V;.rtanA + — 7,, ; -^-tanA 




<?A ib i - 


3k, 


, 2. * ? n ■> > ^y, r A ~ ] /f 

A. cos /\ ' °Ps k: 1 C P< 


I MA 3 — ~ - A 4 fZl: * - ,4- J 

+(4 * 7 6 —^ +b ^ 


?Ps 


i rc L 4 Qh. ~ — A . 5 *y _ t x 

+ C '* 1 ^ ! W y-r‘ anA + i> — ' ■ 1 


/y ( x A 4* b y (< ^ , tan A 


C Px 


*> S - _ r _ _5A_ , , 

2 1 'Uxy'O.x ~ n J^jt 


cos‘A 


<7, 


[5.4.2] 


i 2 > 


r‘M 


d/L 


+[ „,^ tanA 7 ,.^ + ^ f - ^ 


+ (6 4 tan A y (> x7 7,) -~ L + [A 3 y ( y ; , j - tan Ay,. ; y ;> j] ^ 


QPt 


+ ^ * Y ^^r +bi un a *. <? 7 a 
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■ ' 


+ [[^y, y y - tan A \y T 'Mch + j ', a cj 

n 


(' h 4 tan A ^ *y tf A cx + O^t Vj.x-'T' tan Ay, - y ; r ] /* a „ 


+ b 


V>.I 


-4 

X^aa. 


+ i) 5 tan Ay,- y, t ? /f aT ] ]<*>/ 


The derivatives of the displacement function with respect to a general shape parameter 
p l are given as follows: 


<3 _ d dx + <3 dy 

dp, dx dp s dy 8p s 


The global coordinates x and y are given in terms of the local coordinates and the ge- 
ometric shape parameters in equations 5.3.5 and 5.3.6. The derivatives ofx and.v with 
respect to various geometric shape parameters are as follows: 


cy S(>? + 1) 

dAR ~ 4(JaR S , 


[5.4.4] 


dy AR(n + 1) 
dS = 4 (JARS) 


[5.4.5] 


dy 

dtp 


0.0 


[5.4.6] 


<?v 


Ci\ 


- 0.0 


[5.4.7] 
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[3.4.8] 


+ ~r ~ 7 T tan A 


BAR JL ' <3/1 

2(1 + //>)/!/? 5 2 


ax 


dy 

+ T-rrr tan A 


55 2(1 + //>)V7k5 bar 


[5.4.9] 


S^IARS 


d ‘P [(1 +tp)l 2 


[5.4.10] 


; )x Jars (? + i ) 

dt\ 


2 cos‘A 


[5.4.11] 


The derivatives of the ha.'fchord with respect to the shape parameters are: 


db _ (1 + U p - 1 )n -1- tp) cos A 
^ 4(1 + rp)VJRS 


db 5*(l + (tp -l)i/ + tp) cos A 

BAR ~~ ± 

4(1 + tp){S A R) 2 


[5.4.12] 


[5.4.13] 


db Srj cos A 

" S3 ■ " '■■■ ■ »■ 1 " - 

e ‘P (1 + tp) 2 jar's 


[5.4.14] 


db 5(1 + (jp - \)rj + tp) sin A 
” 2(1 + (p)J~S A R 


The aerodynamic coeflicierits, A m , A^, A aA , A n and A„, are functions of the reduced 
frequency. The reduced frequency changes as the shape parameters change therefore 
these terms are functions of the shape parameters and must be included in the formu- 
lation. They are defined as: 
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[5.4.16] 
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[5.4.17] 


[5.4.18] 


[5.4.19] 


[5.4.20] 


[5.4.21] 
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- C(k n ) 


Hr 





dC{k n ) 

dk„ 


k n 


)( 


c, 


/a 

7T 


)]--r^fS.4.22] 


The derivative of Theodorsen's circulation function, as given in equation (3.1.3), is: 


dC(k„) 

dk n 


.0075075/ 


.1005/ 


[ l -~^-] 2 [1 - 


^2 
K J 


[5.4.23] 


The "analytic ' derivative of the eigenvalue with respect to various parameters, namely 
the surface area S, the aspect ratio AR, the taper ratio ip and sweep were compared with 
those obtained using the two previously described methods, namely: (i) the purely finite 
difference method and (ii) a semi-analytic approach in which the desired derivatives were 
obtained using a forward finite difference scheme. The results are showm in Table 5. An 
excellent agreement exists between the various sets of results. However, the values of the 
derivative of the eigenvalue, with respect to sweep, see fourth row of Table 4, obtained 
using the three methods appear to be different from each other. For example, the "ana- 
lytically obtained ' derivative ( case (iii» is about 6.95 percent more than that obtained 
using a purely finite difference approach (case (i)). Similiarly, the value of the same de- 
rivative obtained using a semi-analytic approach (case (ii)), is about 9.36 percent less 
than that obtained using the analytic approach. The reasons for this descrepencies are 
not entirely clear at this stage. However, Figure 28 shows that the analytic method gives 
the best results. 


Figures 18 through 29 show the flutter speed, frequency, and reduced frequency analyt- 
ically based values versus the values which where reanalyzed at different configurations. 
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These are called reanalysis values because .he v.g method was used to solve for the 
Hotter speed a. the perturbed value. These show excellent agreement over a range which 
can be useful to the designer in the intial design phase. 

However, close attention should be paid to an, mode shifting. For instance, in smaller 
aspect ratio configurations the flutter frequency predominant mode changes from the 
third mode to the second mode. This evident in Figure 21. where the duller frequency 
drops for the lower values of the aspect ratio. 
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FI » J 


50.785 a 
49.120 
48.930 
48.200 


60.748 a 
57.449 

57.100 

56.100 


0.9612 b 
0.9193 
0.9147 
0.9000 



8C.850 

78.558 

78.310 

7S.7Q 


119.85 

109.84 

108.83 

109.30 


1.7345 

1.6370 

1.6268 

1.6350 


0.4S12 
0.4153 
0.409 1 
0.4105 


0.1986 

0.1613 

0,1579 

0.1583 


123.63 

118.03 

117.50 

116.0 


2.1397 

2.0429 

2.034 

2.010 


0.3583 

0.3421 

0.3406 

0.3370 



0.1070 


61.993 

61.811 

61.765 

61.307 


53.302 

51.543 

51.394 

51.566 



107.72 


2566 

.2479 

.2471 

.2469 


.04' 
.04' 
.04784 
.04755 


a) These columns represent the finite difference derivatives of the stiffness matrix. 


b) These columns represent the finite difference derivatives of the mass matrix. 


c) Step size used for finite difference calculations 
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Table 5: Comparison of Eigenvalue Derivatives w.r.t Four Parameters 



case (i) 

Finite Difference 

case (ii) 

Semi -Analytic 

case (iii) 
Analytic 

s 

2.4151 E-5 

2 4593 E-5 

2.4571 E-5 

AR 

7.8769 E-6 

8.1694 E-6 

8.1562 E-6 

TP 

2.6719 E-4 

2.6571 E-4 

2.6330 E-4 

SWEEP 

2.8831 E-5 

2.4199 E-5 

2.6761 E-5 
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Figure 19 Rutter Speed vs Surface Area 
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Figure 20 Reduced Frequency vs Surface Area 
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Figure 21 Flutter Frequency vj Aspect Ratio 
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Figure 23 Reduced Frequency v$ Aspect Ratio 
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Figure 24 Flutter Frequency vs Taper Ratio 
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Figure 25 Flutter Speed vs Taper Ratio J 
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Figure 26 Reduced Frequency vs Taper Ratio 
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Figure 27 Rutter Frequency vs Sw*ep Angle 
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Figure 28 Flutter Speed vs Sweep Angle 
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6.0 Conclusions and Suggestions for Future Work 


A method for analyzing the dynamic aeroelastic behavior of a laminated wing has been 
developed. The aerodynamic formulation was taken from Yates' modified finite strip 
method. This was combined with Giles' equtlavent plate model which is capable of ana- 
lyzing cranked wing box structures. Three different approaches were used to obtain the 
sensitivity of the flutter speed, frequency, and reduced frequency. The first was a purley 
numerical approach using finite difference method. The second used the analytic ex- 
pressions for the derivative of the eigenvalue as originally derived by Lancaster [19], and 
a finite difference method to the calculate the derivatives of the aerodynamic, mass and 
stiffness matrices. The third method also used Lancaster's [19] expression for the 
eigenvalue derivative but the derivative of the aerodynamic matrix is computed analyt- 
ically. 

It v as shown that the eigenvalue derivatives for all three cases are in good agreement 
with each other. Also the results for flutter speed and reduced frequency obtained using 
sensitivity based analysis, for a significant range of parameters, are found to be in good 
agreement with those obtained using a complete reanalysis. 
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There is ample room for additional work in the area of Aerodynamic sensitivity analysis. 
This work can be extended to include the analytical derivatives of the mass and stillness 
matrices. This will elcminate any uncertainties in choosing a finite difference step size 
and be of great use to the present optimization codes. 

An improved aerodynamic code could be added to include the effect of chordwise 
bending of the box sta cture. This can oe achieved by using f he aerodynamic coefficients 
developed by Spielberg [26] to modify the present analysis. 

The variation of the curvature of the flutter speed, flutter frequency , and reduced fre- 
quency could be predicted more acurately by a second order derivative. For a general 
complex matrix, efficient calculations of these derivatives are present by Murthy and 
Haftka [19]. 


The aerodynamic sensitivity calculations of a wing with mass stores is needed in the in- 
dustry. The changing mass configurations of a wing with missiles and bombs is a typical 
example for the need of such an analysis capability. 

This work can be extended to include the sensitivities of the flutter speed, flutter fre- 
quency, and reduced frequency to a simultaneous change in several shape parameters. 
For example, the sensitivity of flutter speed with the changes in sweep and aspect ratio 
would be useful. The structural code used in this analysis is extremely versatile and can 
be extended with the present aerodynamic formulation to analyze cranked wing boxes. 
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